April  2014 


CARRIER  ET  AL. 


1509 


Impact  of  Assimilating  Ocean  Velocity  Observations  Inferred  from  Lagrangian 
Drifter  Data  Using  the  NCOM-4DVAR* 

Matthew  J.  Carrier,  Hans  Ngodock,  Scott  Smith,  and  Gregg  Jacobs 

Naval  Research  Laboratory,  Stennis  Space  Center,  Mississippi 

Philip  Muscarella 

American  Society  for  Engineering  Education,  Washington,  D.  C. 

Tamay  Ozgokmen  and  Brian  Haus 

University  of  Miami,  Miami,  Florida 

Bruce  Lipphardt 

University  of  Delaware,  Newark,  Delaware 

(Manuscript  received  23  July  2013,  in  final  form  6  November  2013) 

ABSTRACT 

Eulerian  velocity  fields  are  derived  from  300  drifters  released  in  the  Gulf  of  Mexico  by  The  Consortium  for 
Advanced  Research  on  Transport  of  Hydrocarbon  in  the  Environment  (CARTHE)  during  the  summer  2012 
Grand  Lagrangian  Deployment  (GLAD)  experiment.  These  data  are  directly  assimilated  into  the  Navy 
Coastal  Ocean  Model  (NCOM)  four-dimensional  variational  data  assimilation  (4DVAR)  analysis  system  in 
a  series  of  experiments  to  investigate  their  impact  on  the  model  circulation.  The  NCOM-4DVAR  is  a  newly 
developed  tool  for  data  analysis,  formulated  for  weak-constraint  data  assimilation  based  on  the  indirect 
representer  method.  The  assimilation  experiments  take  advantage  of  this  velocity  data  along  with  other 
available  data  sources  from  in  situ  and  satellite  measurements  of  surface  and  subsurface  temperature  and 
salinity.  Three  different  experiments  are  done:  (i)  A  nonassimilative  NCOM  free  run,  (ii)  an  assimilative 
NCOM  run  that  utilizes  temperature  and  salinity  observations,  and  (iii)  an  assimilative  NCOM  run  that  uses 
temperature  and  salinity  observations  as  well  as  the  GLAD  velocity  observations.  The  resulting  analyses  and 
subsequent  forecasts  are  compared  to  assimilated  and  future  GLAD  velocity  and  temperature/salinity  ob¬ 
servations  to  determine  the  performance  of  each  experiment  and  the  impact  of  the  GLAD  data  on 
the  analysis  and  the  forecast.  It  is  shown  that  the  NCOM-4DVAR  is  able  to  fit  the  observations  not  only  in  the 
analysis  step,  but  also  in  the  subsequent  forecast.  It  is  also  found  that  the  GLAD  velocity  data  greatly  im¬ 
proves  the  characterization  of  the  circulation,  with  the  forecast  showing  a  better  fit  to  future  GLAD  obser¬ 
vations  than  those  experiments  without  the  velocity  data  included. 


1.  Introduction 

Modern  ocean  forecasting  using  numerical  models  has 
advanced  significantly  in  the  past  20  years.  Models,  such 
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as  the  Princeton  Ocean  Model  (POM;  Blumberg  and 
Mellor  1987),  the  Regional  Ocean  Modeling  System 
(ROMS;  Shchepetkin  and  McWilliams  2003,  2005; 
Marchesiello  et  al.  2001),  the  Navy  Coastal  Ocean 
Model  (NCOM;  Martin  2000;  Barron  et  al.  2006),  and 
the  Hybrid  Coordinate  Ocean  Model  (HYCOM;  Bleck 
2002)  are  capable  of  producing  very  accurate  ocean 
simulations  for  regional  and  global  applications  at  reso¬ 
lutions  ranging  from  tens  of  kilometers  to  a  few  hundred 
meters.  Despite  recent  advancement,  numerical  model¬ 
ing  continues  to  suffer  from  inherent  errors  associated 
with  unresolved  processes,  inaccurate  parameterizations, 


©  2014  American  Meteorological  Society 


Report  Documentation  Page 


Form  Approved 
OMB  No.  0704-0188 


Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 
VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number. 


1.  REPORT  DATE 

APR  2014 


2.  REPORT  TYPE 


4.  TITLE  AND  SUBTITLE 

Impact  of  Assimilating  Ocean  Velocity  Observations  Inferred  from 
Lagrangian  Drifter  Data  Using  the  NCOM-4DVAR 


6.  AUTHOR(S) 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

Naval  Research  Laboratory, Bldg.  1009,  Balch  Blvd,Stennis  Space 
Center, MS, 39529 

9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 


3.  DATES  COVERED 

00-00-2014  to  00-00-2014 

5a.  CONTRACT  NUMBER 

5b.  GRANT  NUMBER 

5c.  PROGRAM  ELEMENT  NUMBER 

5d.  PROJECT  NUMBER 

5e.  TASK  NUMBER 

5f.  WORK  UNIT  NUMBER 

8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 


10.  SPONSOR/MONITOR’S  ACRONYM(S) 

11.  SPONSOR/MONITOR’S  REPORT 
NUMBER(S) 


12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  unlimited 

13.  SUPPLEMENTARY  NOTES 

14.  ABSTRACT 

Eulerian  velocity  fields  are  derived  from  300  drifters  released  in  the  Gulf  of  Mexico  by  The  Consortium  for 
Advanced  Research  on  Transport  of  Hydrocarbon  in  the  Environment  (CARTHE)  during  the  summer 
2012  Grand  Lagrangian  Deployment  (GLAD)  experiment.  These  data  are  directly  assimilated  into  the 
Navy  Coastal  Ocean  Model  (NCOM)  four-dimensional  variational  data  assimilation  (4DVAR)  analysis 
system  in  a  series  of  experiments  to  investigate  their  impact  on  the  model  circulation.  The  NCOM-4DVAR 
is  a  newly  developed  tool  for  data  analysis,  formulated  for  weak-constraint  data  assimilation  based  on  the 
indirect  representer  method.  The  assimilation  experiments  take  advantage  of  this  velocity  data  along  with 
other  available  data  sources  from  in  situ  and  satellite  measurements  of  surface  and  subsurface  temperature 
and  salinity.  Three  different  experiments  are  done:  (i)  A  nonassimilative  NCOM  free  run,  (ii)  an 
assimilative  NCOM  run  that  utilizes  temperature  and  salinity  observations,  and  (iii)  an  assimilative  NCOM 
run  that  uses  temperature  and  salinity  observations  as  well  as  the  GLAD  velocity  observations.  The 
resulting  analyses  and  subsequent  forecasts  are  compared  to  assimilated  and  future  GLAD  velocity  and 
temperature/salinity  observations  to  determine  the  performance  of  each  experiment  and  the  impact  of  the 
GLAD  data  on  the  analysis  and  the  forecast.  It  is  shown  that  the  NCOM-4DVARis  able  to  fit  the 
observations  not  only  in  the  analysis  step,  but  also  in  the  subsequent  forecast.  It  is  also  found  that  the 
GLAD  velocity  data  greatly  improves  the  characterization  of  the  circulation,  with  the  forecast  showing  a 
better  fit  to  future  GLAD  observations  than  those  experiments  without  the  velocity  data  included. 


15.  SUBJECT  TERMS 


16.  SECURITY  CLASSIFICATION  OF: 


a.  REPORT 

unclassified 


b.  ABSTRACT 

unclassified 


c.  THIS  PAGE 

unclassified 


17.  LIMITATION  OF 

18.  NUMBER 

ABSTRACT 

OF  PAGES 

Same  as 

16 

Report  (SAR) 

19a.  NAME  OF 
RESPONSIBLE  PERSON 


Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


1510 


MONTHLY  WEATHER  REVIEW 


Volume  142 


errors  in  specified  surface  and/or  boundary  conditions, 
and  errors  in  the  initial  conditions.  As  such,  if  a  model 
simulation  were  to  continue  uncorrected,  the  forecast 
trajectory  will  inevitably  veer  away  from  reality  as  the 
effects  of  these  errors  accumulate. 

Because  of  the  impact  of  model  error,  frequent  updating 
of  the  model  trajectory  from  available  observations  is 
a  key  component  of  any  ocean  or  atmospheric  fore¬ 
casting  system.  Data  assimilation  (DA)  provides  the 
means  to  produce  suitable  state  estimation  from  which 
the  model  can  be  periodically  adjusted  to  move  the 
forecast  trajectory  closer  to  reality.  For  ocean  fore¬ 
casting,  this  need  is  important  when  attempting  to  ac¬ 
curately  simulate  ocean  currents  at  regional  scales.  On 
large  scales,  existing  ocean  models  do  an  adequate  job  of 
resolving  the  oceanic  flow  from  the  placement  of  the 
Gulf  Stream  to  the  Kuroshio.  However,  at  regional  and 
mesoscales,  where  the  prediction  of  individual  tracers  or 
drifters  are  important  for  search  and  rescue  operations 
or  hydrocarbon/chemical  spill  simulations,  the  ocean 
models  lack  sufficient  accuracy;  therefore,  it  is  important 
to  be  able  to  properly  constrain  the  model  simulation,  at 
the  prescribed  resolution,  with  actual  observations  of 
ocean  current  velocities. 

There  has  been  some  work  in  the  area  of  assimilating 
ocean  velocity  observations  to  improve  the  character¬ 
ization  of  model  currents.  These  have  taken  the  form  of 
simple  nudging  methods  (Fan  et  al.  2004)  to  the  more 
complex  variational  methods  (Taillandier  et  al.  2006; 
Nilsson  et  al.  2012).  Ocean  velocity  DA  has  also  seen 
variety  regarding  the  form  of  the  observation  to  be  used, 
specifically  Lagrangian  (i.e.,  following  the  flow)  or 
Eulerian  (stationary).  Traditional  Eulerian  measurements 
consist  of  point  measurements  from  acoustic  Doppler 
current  profilers  (ADCP)  or  high-frequency  (HF)  radar 
surface  current  measurements  that  measure  the  speed 
and  direction  of  the  ocean  currents  at  a  fixed  location  in 
space.  Lagrangian  data  are  collected  from  data-gathering 
devices  on  board  any  passive  tracer,  such  as  drifters  or 
surface  floats.  Eulerian  observations  can  be  assimilated 
directly  into  the  ocean  model  as  the  form  of  the  data 
matches  that  of  the  model  variable.  Lagrangian  data,  on 
the  other  hand,  must  be  handled  differently  as  the  data 
here  consists  of  position — time  measurements,  not  ocean 
velocity  directly.  There  have  been  numerous  attempts  to 
assimilate  Lagrangian  data,  either  by  deriving  Eulerian 
velocity  from  the  dataset  by  calculating  the  change  in 
drifter  position  over  some  time  scale  (Hernandez  et  al. 
1995;  Ishikawa  et  al.  1996),  or  by  assimilating  the  La¬ 
grangian  data  directly  by  evolving  a  series  of  tracers 
within  the  model  simulation  then  minimizing  the  distance 
between  the  observed  and  model  tracers  via  some  DA 
scheme  (Molcard  et  al.  2003;  Ozgokmen  et  al.  2003; 


Taillandier  et  al.  2006;  Nilsson  et  al.  2012),  or  by  state 
augmentation  via  a  Kalman  filter  scheme  (Ide  et  al.  2002; 
Kuznetsov  et  al.  2003;  Salman  et  al.  2006).  Molcard  et  al. 
(2005)  point  out  that  the  former  method,  dubbed 
“pseudo-Lagrangian,”  can  be  problematic  if  the  data 
sampling  period  exceeds  the  Lagrangian  time  scale, 
which  is  usually  on  the  order  of  1-3  days  for  ocean 
surface  velocity  measurements  and  7-15  days  for  the 
ocean  interior. 

For  this  work,  deriving  Eulerian  velocity  measurements 
from  the  drifter  positions  from  the  Grand  Lagrangian 
Deployment  (GLAD)  experiment  is  not  problematic,  as 
the  sampling  time  is  every  5  min  (Pojc  et  al.  2013, 
manuscript  submitted  to  Nature );  therefore,  derivation 
of  reasonably  accurate  Eulerian  velocities  is  possible. 
This  effort  aims  to  assess  the  impact  of  Eulerian  velocity 
data,  derived  from  the  GLAD  drifter  locations,  on 
a  simulation  of  the  Gulf  of  Mexico  using  NCOM.  These 
observations  are  assimilated  via  the  four-dimensional 
variational  data  assimilation  (4DVAR)  representer 
method  (Bennett  1992,  2002;  Chua  and  Bennett  2001; 
Kurapov  et  al.  2009,  2011)  employed  by  the  NCOM- 
4DVAR  (Ngodock  and  Carrier  2013).  These  observations, 
along  with  available  remote  and  in  situ  observations  of 
temperature  and  salinity  are  used  to  initialize  a  series  of 
4-day  forecasts.  These  analyses  and  forecasts  are  exam¬ 
ined  by  comparison  to  future  GLAD -derived  Eulerian 
velocity  observations  and  traditional  temperature/salinity 
observations  in  an  attempt  to  evaluate  the  impact  of  the 
assimilated  velocities  on  the  analyses  and  subsequent 
forecasts. 

This  paper  is  organized  in  the  following  manner:  sec¬ 
tion  2  provides  a  discussion  of  the  forecast  model  and 
the  accompanying  assimilation  system;  section  3  pro¬ 
vides  an  overview  of  the  GLAD  experiment  as  well  as 
the  data  processing  methodology  used  to  derive  Eulerian 
values  from  the  drifter  position  data;  section  4  presents 
the  design  of  the  experiments  used  to  evaluate  the  impact 
of  the  GLAD  data  on  the  ocean  analysis  and  forecast; 
section  5  provides  a  discussion  of  the  experiment  results 
followed  by  a  general  summary  and  a  discussion  of  fu¬ 
ture  work  in  section  6. 

2.  Forecast  model  and  analysis  system 

a.  The  Navy  Coastal  Ocean  Model  (NCOM) 

NCOM  is  selected  as  the  ocean  forecast  model  for  this 
work  and  is  capable  of  producing  forecasts  of  tempera¬ 
ture,  salinity,  ocean  currents,  and  sea  surface  height  for 
regional  and  global  applications.  It  is  a  primitive  equa¬ 
tion  model  with  a  free-surface  and  a  generalized  vertical 
coordinate  that  can  be  configured  with  terrain-following 
free-sigma  or  fixed-sigma,  or  constant  z-level  surfaces  in 
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a  number  of  combinations  (Barron  et  al.  2006).  The 
model  employs  the  Mellor-Yamada  level-2.5  turbu¬ 
lence  closure  parameterization  (Mellor  and  Yamada 
1982)  for  vertical  diffusion  and  the  Smagorinsky  scheme 
(Smagorinsky  1963)  for  horizontal  diffusion. 

The  model  domain  for  this  experiment  extends  from 
18°-31°N  to  79°-98°W  using  a  spherical  coordinate 
projection  at  a  horizontal  resolution  of  6  km.  The  model 
has  50  layers  in  the  vertical,  with  25  free-sigma  levels 
extending  to  a  depth  of  116  m  with  constant  z-levels 
extending  down  to  a  maximum  of  5500  m  with  the  depth 
of  the  first  subsurface  layer  at  0.5  m.  The  model  resolu¬ 
tion  is  coarse  compared  to  other  simulations  of  the  Gulf 
of  Mexico,  but  was  necessary  because  of  the  computa¬ 
tional  expense  of  the  NCOM-4DVAR,  which  can  be 
many  times  more  expensive  than  the  forward  model 
(Ngodock  and  Carrier  2013).  Lateral  boundary  condi¬ 
tions  are  provided  by  the  global  NCOM  model  at  1/8° 
resolution  (every  3  h)  and  surface  atmospheric  forcing, 
such  as  wind  stress,  atmospheric  pressure,  and  surface 
heat  flux  is  provided  by  the  0.5°  Navy  Operational 
Global  Atmospheric  Prediction  System  (NOGAPS) 
model  every  3h  (Rosmond  et  al.  2002);  river  forcing  is 
provided  at  all  river  in-flow  locations  in  the  Gulf  of 
Mexico  domain. 

b.  The  NCOM  four-dimensional  variational 
assimilation  system  ( NCOM-4D  VAR ) 

The  DA  system  selected  for  this  work  is  the  NCOM- 
4DVAR  developed  using  the  dynamical  core  of  the 
forward  NCOM  model.  The  NCOM-4DVAR  is  a  vari¬ 
ational  assimilation  system  based  on  the  indirect  rep¬ 
resenter  method  as  described  by  Bennett  (1992,  2002) 
and  Chua  and  Bennett  (2001).  This  system  has  been 
described  in  detail  by  Ngodock  and  Carrier  (2013),  and 
a  full  derivation  of  the  representer  method  can  be  found 
in  Chua  and  Bennett  (2001),  therefore,  only  an  overview 
is  provided  here.  The  representer  method  aims  to  find  an 
optimal  analysis  solution  as  the  linear  combination  of 
a  first  guess  (i.e.,  prior  model  solution)  and  a  finite 
number  of  representer  functions,  one  per  datum 

M 

u(x,t)  =  up(x,t)+  X  Pmrm(x’t),  (1) 

m— 1 

where  u(x ,  t )  is  the  optimal  analysis  solution,  up(x,  t)  is 
the  prior  forecast,  rm(x ,  t )  is  the  representer  function  for 
the  rath  observation,  and  j3m  is  the  rath  representer  co¬ 
efficient.  The  representer  coefficients  can  be  found  by 
solving  the  linear  system: 

(R  +  0)/3  =  y— HxL  (2) 


where  O  is  the  observation  error  covariance,  y  is  the 
observation  vector,  H  is  the  linear  observation  operator 
that  maps  the  model  fields  to  the  observation  locations, 

is  the  model  vector,  and  R  is  the  representer  matrix 
and  is  equivalent  to  HMBMTHT  (M  is  the  tangent  linear 
model,  or  TLM;  MT  is  the  adjoint  of  NCOM;  B  is  the 
model  error  covariance;  and  T  denotes  the  linear 
transposition).  Since  the  matrix  R  +  O  is  symmetric  and 
positive  definite,  (2)  can  be  solved  for  /3  iteratively  using 
a  linear  solver,  such  as  the  conjugate  gradient  method. 
From  (2)  it  is  clear  that  the  for  each  representer  can 
be  found  by  integrating  the  adjoint  and  TLM  over  some 
number  of  minimization  steps  until  convergence.  Once 
found,  Pm  is  acted  upon  by  (1),  involving  one  final 
sweep  through  the  adjoint  and  TLM  to  find  the  optimal 
correction. 

In  the  NCOM-4DVAR,  j3m  is  found  with  a  precondi¬ 
tioned  conjugate  gradient  solver.  The  preconditioner 
here  follows  from  Courtier  (1997)  to  introduce  a  change 
of  variable  in  the  minimization  step  described  in  (2), 
where  /3  is  redefined  as  u  =  VO/3  so  that  (2)  can  now  be 
expressed  as 

(VOrTRVO::T+l)u=  v/0=I(y-Hx/).  (3) 

This  transformation  ensures  that  there  is  a  lower  bound 
of  1  for  the  eigenvalues,  which  ensures  that  the  condition 
number  will  remain  reasonably  small. 

The  background  and  model  error  covariance  in 
NCOM-4DVAR  follow  the  work  of  Weaver  and  Courtier 
(2001)  and  Carrier  and  Ngodock  (2010)  where  the  error 
covariance  is  univariate.  This  is  deemed  acceptable  as 
the  application  of  the  tangent  linear  and  adjoint  models 
in  the  minimization  and  final  sweep  provide  multivariate 
balance  constraints  through  the  linearized  dynamics.  It 
has  been  shown  (Yu  et  al.  2012)  that  omitting  linear 
balance  constraints  does  not  lead  to  a  significant  deg¬ 
radation  of  the  final  solution  in  terms  of  the  fit  to  ob¬ 
servations.  The  univariate  error  covariance  can  be 
further  decomposed  into  a  correlation  matrix  and  the 
associated  error  variance  such  that 

B  =  2C2,  (4) 

where  2  is  a  diagonal  matrix  of  the  error  standard  de¬ 
viation  and  C  is  a  symmetric  matrix  of  error  correlations. 
In  NCOM-4DVAR  the  error  standard  deviations  of  the 
background  are  used  at  the  initialization  of  the  tangent 
linear  model  only,  whereas  the  model  error  (also  con¬ 
tained  in  the  matrix  2)  is  used  when  the  adjoint  forces 
the  tangent  linear  model  during  integration  (i.e.,  as  the 
tangent  linear  model  integrates  forward  in  time).  This 
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allows  the  weak  constraint  method  to  correct  for  the 
initial  condition  error  while  also  adjusting  the  forward 
model  trajectory  based  on  the  specification  of  the  model 
error.  The  error  correlation,  for  both  the  model  and  the 
background  errors  are  not  directly  calculated  and  stored 
in  NCOM-4DVAR;  rather,  the  effect  of  the  correlation 
matrix  acting  on  an  input  vector  is  modeled  by  the  so¬ 
lution  of  a  diffusion  equation  following  the  work  of 
Weaver  and  Courtier  (2001).  For  a  full  explanation  of 
this  method,  we  refer  the  reader  to  Weaver  and  Courtier 
(2001)  or  Yaremchuk  et  al.  (2013);  for  a  description  of 
the  implementation  of  this  method  in  NCOM-4DVAR, 
we  refer  the  reader  to  Carrier  and  Ngodock  (2010)  or 
Ngodock  (2005). 

3.  Data  and  processing 

With  a  Gulf  of  Mexico  Research  Initiative  (GoMRI) 
award,  the  Consortium  for  Advanced  Research  on  Trans¬ 
port  of  Hydrocarbons  in  the  Environment  (CARTHE) 
conducted  an  unprecedented  large-scale  deployment  of 
300  custom-built  drifters  in  the  northern  Gulf  of  Mexico 
in  the  summer  of  2012.  Equipped  with  GPS  positioning, 
the  drifters  are  capable  of  reporting  their  positions  every 
5  min,  allowing  for  excellent  temporal  resolution  and  for 
accurate  estimates  of  the  Eulerian  velocity  along  the 
drifter  track.  The  drifters  are  drogued  to  a  depth  of  1  m 
(Pojc  et  al.  2013,  manuscript  submitted  to  Nature). 

Although  the  NCOM-4DVAR  system  can  be  adapted 
to  assimilate  drifter  positions  because  of  the  available 
adjoint  of  the  momentum  equation,  it  is  decided  for 
the  sake  of  convenience  to  use  the  drifter  data  to  esti¬ 
mate  the  Eulerian  velocities,  which  could  be  more  easily 
assimilated  into  the  model  using  the  existing  assimila¬ 
tion  architecture  without  the  need  to  update  or  alter  the 
existing  forward,  tangent  linear,  and  adjoint  codes.  The 
basic  principle  behind  the  derivation  of  Eulerian  ve¬ 
locity  from  a  Lagrangian  data  source  involves  the  cal¬ 
culation  of  the  distance  traveled  by  the  drifter  over  some 
time  scale.  In  practice,  however,  the  calculation  is  more 
complicated  due  to  noisy  position  values  and  non¬ 
physical  accelerations,  such  as  a  drifter  being  removed 
from  the  water  by  the  crew  of  a  passing  ship.  The  GLAD 
drifter  observations  are  highly  nonstationary,  and  there¬ 
fore,  can  simultaneously  represent  multiple  scales  at 
once.  This  makes  filtering  the  observations  very  difficult, 
because  filtering  the  larger  scales  of  motion  will  overly 
smooth  the  small-scale  velocity  values.  Therefore, 
a  multistep  process  is  used  to  derive  Eulerian  velocities 
from  the  GLAD  drifter  position  data.  First,  the  in¬ 
stantaneous  velocity  and  its  associated  acceleration 
vector  (i.e.,  change  in  instantaneous  velocity  from  one 
velocity  location  to  the  next)  are  computed.  The  algorithm 


cancels  any  data  point  where  the  velocity  magnitude 
exceeds  the  threshold  maximum  of  3  m  s-1  or  if  the  data 
rotation  exceeds  the  maximum  threshold  of  2tt  (3h)_1. 
No  subsequent  data  from  a  particular  drifter  are  in¬ 
cluded  in  the  processing  if  the  drifter  travels  more  than 
80  km  in  12  h,  or  if  the  drifter  travels  less  than  100  m  in 
12  h.  The  remaining  velocities  are  sampled  at  the 
smallest  possible  constant  interval  and  are  filtered  to 
remove  high-frequency  noise  (cutoff  of  2  X  Nyquist, 
roughly  1/10  min).  The  velocity  values  are  then  com¬ 
puted  at  each  remaining  position  and  the  data  are 
stored. 

In  addition  to  the  derived  Eulerian  velocity  data, 
other  available  sources  of  data  are  used  in  this  study. 
Both  remotely  sensed  and  in  situ  ocean  observation  data 
are  assimilated  from  (Geostationary  Operational  Envi¬ 
ronmental  Satellite)  GOES-East  sea  surface  tempera¬ 
tures  (SSTs),  Argo  profiling  floats  (Roemmich  et  al. 
2001),  expendable  bathythermographs  (XBT),  and 
drifting  buoys.  These  data  are  gathered  and  quality 
controlled  using  the  operational  data  preparation  utility 
from  the  Navy  Coupled  Ocean  Data  Assimilation 
(NCODA)  system  (Cummings  2005).  These  data  are 
collected,  preprocessed  and  used  within  the  NCOM- 
4DVAR  assimilation  window  at  their  respective  obser¬ 
vation  times.  It  should  be  noted  that  no  altimeter  data 
are  assimilated  in  this  work.  This  is  done  in  order  to 
properly  evaluate  the  impact  of  the  velocity  observations 
in  the  assimilation.  Since  altimeter  observations  help 
constrain  the  mesoscale  features,  it  can  be  used  to  correct 
ocean  velocity  as  well.  In  this  present  work,  corrections 
made  to  the  mesoscale  features  and  the  ocean  velocity 
can  be  directly  attributed  to  the  assimilation  of  the  drifter 
observations  alone. 

4.  Experiment  design 

The  forecast  and  analysis  experiments  are  conducted 
to  ascertain  the  impact  of  the  GLAD  data  on  the  NCOM 
forecast  and  to  determine  if  their  inclusion  within  the 
data  stream  improves  the  forecast  of  near-surface  ocean 
currents.  The  experiments  shown  in  this  work  cover  the 
time  period  between  1  August  and  30  September  2012 
and  are  located  in  the  region  where  the  GLAD  drifters 
are  deployed  during  this  time  frame.  During  this  time, 
the  loop  current  has  already  shed  an  eddy  that  is  located, 
for  the  majority  of  the  experiment,  near  the  central  Gulf 
of  Mexico.  Figure  1  (left  panel)  shows  a  composite  im¬ 
age  of  sea  surface  height  (SSH)  observations  taken  from 
the  available  altimeters  between  20  and  30  September 
2012  and  processed  by  the  NCODA  data  preparatory 
suite.  The  estimated  heights  are  shown  here,  not  the 
anomalies,  for  a  more  direct  comparison  to  the  model 
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Fig.  1.  (left)  Composite  of  observed  sea  surface  height  (SSH)  values  from  20  to  30  Sep  2012  (value  of  the  observation 
indicated  by  color  bar  in  meters),  (right)  NCOM  model  SSH  solution,  valid  30  Sep  2012. 


SSH  field,  shown  in  Fig.  1  (right  panel).  This  model  SSH 
field  is  valid  on  30  September  2012  (taken  from  a  free- 
run  NCOM  solution).  To  estimate  the  sea  surface  height 
from  the  observation  data,  a  long-term  mean  SSH  field 
from  HYCOM  has  been  interpolated  to  the  observation 
locations  and  added  to  the  anomalies.  Figure  1  indicates 
that  the  loop  current  has  shed  an  eddy,  now  located  near 
the  central  Gulf  of  Mexico.  The  NCOM  solution  in 
Fig.  1  shows  that  the  forecast  model  has  captured  the 
general  location  of  the  eddy  and  its  separation  from  the 
loop  current,  although  the  shape  and  positioning  does 
not  exactly  match  the  altimeter  data  (its  position  being 
offset  to  the  west  from  the  observations). 

The  location  of  the  GLAD  drifter  observations  covers 
much  of  the  central  Gulf  of  Mexico  throughout  the  ex¬ 
periment  time  frame.  Figure  2  shows  the  location  of  the 
processed  GLAD  drifter  velocity  data  from  1  August  to 
30  September  2012  (positions  plotted  daily).  It  can  be 
seen  in  Fig.  2  that  the  GLAD  drifter  velocity  data  has 
a  high  concentration  of  observations  near  the  northern 
Gulf  Coast,  but  still  with  good  coverage  farther  south 
near  the  loop  current  and  the  region  of  the  loop  current 
eddy. 

a.  Analysis  procedure 

For  this  study,  the  observations  are  assimilated  in 
a  weak-constraint  mode  in  order  to  have  full  control 
over  the  correction  of  the  model  trajectory  and  to  at¬ 
tempt  to  account  for  sources  of  model  error  that  we 
believe  have  the  largest  contribution  to  the  forecast  er¬ 
ror  in  surface  current  velocity.  The  initial  condition  er¬ 
ror  for  each  of  the  following  experiments  is  set  as  2.0°C 
for  temperature,  0.5  practical  salinity  unit  (psu)  for  sa¬ 
linity,  0.2 ms-1  for  velocity,  and  0.1m  for  SSH.  These 
errors  have  been  set  by  examining  the  innovation  values 
of  the  96-h  forecasts  from  the  free-running  NCOM 
model  as  compared  to  available  observations.  The  96-h 
forecast  is  used  to  estimate  the  initial  condition  error 


due  to  the  update  cycle  employed  in  this  work  (see 
section  4b).  These  values  are  uniform  across  the  ocean 
model  domain,  mainly  to  aid  in  the  convergence  rate  of 
the  assimilation  system,  but  also  the  values  set  here  are 
representative  of  the  average  error  across  the  domain. 
These  values  are  reduced  at  depth  as  exp(z/zo),  where 
zo  —  1000  m  and  z  is  0  at  the  surface  and  negative  at 
depth  (Yu  et  al.  2012).  This  is  deemed  acceptable  as  it  is 
not  expected  that  the  ocean  model  will  be  largely  in 
error  at  depth,  especially  considering  that  we  are  mostly 
interested  in  the  accuracy  of  surface  currents.  The  model 
errors  are  assumed  to  be  due  to  errors  in  the  specified 
atmospheric  surface  forcing  for  this  study.  This  is 
deemed  reasonable  as  ocean  currents  at  the  surface  are 
strongly  affected  by  surface  wind  stress.  As  such,  the 
errors  are  assumed  to  be  10%  of  the  average  value  of 
each  surface  forcing  field  through  the  assimilation  win¬ 
dow,  which  corresponds  to  roughly  40Wm”2  for  the 
combined  solar  and  temperature  fluxes,  0.25  Pa  for 
the  surface  wind  stress,  and  4.0  X  10-7  psu  m  s-1  for  the 
surface  salinity  flux.  These  error  values  are  then  con¬ 
verted  to  the  units  of  the  ocean  model  variables  (i.e., 


Fig.  2.  Location  of  each  GLAD  drifter  velocity  observation  from 
1  Aug  to  30  Sep  2012  (observations  plotted  at  daily  intervals). 
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heat  flux  to  ocean  temperature,  wind  stress  to  ocean 
velocity)  and  applied  directly  to  the  ocean  variables  of 
temperature,  salinity,  and  velocity.  The  scale  of  the 
initial  and  model  errors  is  taken  to  be  roughly  equivalent 
to  the  Rossby  radius  of  deformation,  which  for  the  Gulf 
of  Mexico  ranges  between  30  and  40  km  and  is  fixed  in 
time.  Although  the  background  error  scales  are  static  in 
time,  the  time-evolving,  flow-dependent  scales  are  ap¬ 
plied  via  the  dynamics  provided  by  the  tangent  linear 
and  adjoint  of  NCOM. 

Each  of  the  observations  is  also  assigned  errors.  In 
operational  data  assimilation,  the  observation  error  is 
a  combination  of  the  estimated  instrument  error  and  the 
assumed  representative  error,  with  the  latter  set  to  ac¬ 
count  for  features  that  may  be  present  in  the  dataset  that 
the  model  does  not  resolve.  For  this  experiment  each 
observation  is  assigned  a  static  error  that  does  not  vary 
in  space  or  time,  again  to  aid  in  the  convergence  of  the 
assimilation  system.  The  temperature  and  salinity  ob¬ 
servation  errors  are  representative  of  the  values  pro¬ 
duced  by  the  operational  data  processing  system 
employed  by  the  U.S.  Navy.  The  errors  are  set  as  0.2°C 
for  temperature,  0.1  psu  for  salinity,  and  0.02 ms-1  for 
velocity.  The  velocity  error  is  found  as  a  consequence 
of  the  GLAD  data  processing.  The  dataset  used  to 
compute  the  velocities  for  this  work  was  the  15-min 
interval  dataset  derived  from  the  GLAD  position  data. 
The  drifter  position  error  is  estimated  at  10  m.  The 
corresponding  velocity  error  is  taken  as  the  sum  of  the 
two  position  data  points  (—20  m)  divided  by  the  time 
interval  between  them  (900  s).  This  yields  a  velocity 
error  of  about  0.02 ms-1. 

b.  Selected  experiments 

Each  experiment  proceeds  as  a  series  of  4-day  win¬ 
dows  from  1  August  to  30  September  2012.  At  the  end  of 
each  4-day  assimilation  window,  the  forecast  model  is 
run  from  the  updated  initial  condition  to  provide  the 
background  for  the  next  4-day  assimilation  period.  A 
free-run  forecast  (i.e.,  no  assimilation)  is  also  run  in 
order  to  compare  the  assimilation  analysis  and  forecast 
to  a  control  to  evaluate  the  impact  of  the  assimilated 
observations.  In  a  series  of  preliminary  experiments  (not 
shown),  an  attempt  was  made  to  evaluate  the  impact  of 
the  temporal  frequency  at  which  the  observations  are 
sampled  for  the  assimilation.  Three  individual  tests  were 
done,  each  consisting  of  one  4-day  assimilation  cycle  and 
subsequent  4-day  forecast  in  which  only  the  data  sam¬ 
pling  is  varied  while  ingesting  available  observations 
every  6  h,  3  h,  and  hourly  to  determine  which  is  the  best 
setting  for  this  study.  It  was  determined  that  sampling 
observations  every  hour  provides  the  best  quality  in  terms 
of  analysis  fit  to  the  observations  as  well  as  subsequent 


improvement  to  the  forecast;  the  increased  computa¬ 
tional  burden  over  less  frequent  observation  input  is 
minimal.  It  is  interesting  to  note  that  this  test  illuminates 
the  need  for  a  larger  number  of  velocity  measurements 
in  the  assimilation  step  to  ensure  a  good  analysis  and 
forecast  fit  to  the  observations.  This  is  not  surprising, 
however,  as  the  ocean  velocity  is  highly  variable  in  time 
and  space;  therefore,  more  data  are  needed  to  properly 
constrain  the  model.  All  subsequent  experiments  shown, 
therefore,  employ  the  hourly  sampling  strategy. 

Two  primary  experiments  are  carried  out  to  evaluate 
the  impact  of  the  GLAD  velocity  observations  on 
the  NCOM-4DVAR  analysis  and  subsequent  NCOM 
forecast.  These  experiments  are  (i)  a  cycling  analysis/ 
forecast  run  from  1  August  to  30  September  2012  that 
uses  the  NCOM-4DVAR  to  assimilate  temperature  and 
salinity  observations  only  (hereafter  referred  to  as  TS) 
and  (ii)  a  cycling  analysis/forecast  run  from  1  August  to 
30  September  where  temperature  and  salinity  observa¬ 
tions  are  assimilated  along  with  the  GLAD  velocity 
observations  by  the  NCOM-4DVAR  analysis  system 
(hereafter  referred  to  as  ALL).  The  analyses  of  each  of 
these  experiments  are  compared  to  the  assimilated  data 
to  evaluate  the  performance  of  the  assimilation  in  terms 
of  the  fit  to  the  observations.  The  subsequent  forecast 
from  each  analysis  is  compared  to  unassimilated  data  in 
the  next  analysis  window  (data  that  have  yet  to  be  as¬ 
similated,  i.e.,  examining  the  innovations  prior  to  as¬ 
similation).  Theoretically,  as  data  are  assimilated  by  the 
system,  the  subsequent  forecasts  should  be  improved 
and,  therefore,  the  size  of  the  innovations  in  each  sub¬ 
sequent  assimilation  window  will  be  smaller;  this  is  ex¬ 
amined  in  this  study. 

5.  Results 

As  a  first-order  evaluation  of  the  impact  of  velocity 
data  on  the  assimilation,  the  analysis  fit  to  the  assimi¬ 
lated  observations  is  examined.  Typically,  a  variational 
assimilation  system  should  fit  each  of  the  assimilated 
observations  to  within  one  standard  deviation  of  the 
observation  error.  To  determine  whether  this  criteria  is 
met  for  this  work,  a  normalized  error  metric  is  in¬ 
troduced  as 


=  1  y  \sm-Hm*a\ 

i  <r 


where  ym  is  the  rath  observation,  Xa  is  the  analysis  vector 
in  model  space  (Hm  maps  the  model  analysis  to  the 
observation  location),  am  is  the  error  standard  deviation 
of  the  rath  observation,  and  M  is  the  total  number  of 
observations.  By  (5),  if  the  analysis  residual  is  within  the 
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Fig.  3.  The  /fit  metric  values  for  the  NCOM  free-run  (FR)  model  solution  (solid  line)  and  TS 
analysis  solution  (dashed  line)  measured  against  (top)  assimilated  temperature  observations, 
(middle)  salinity  observations,  and  (bottom)  GLAD  velocity  observations.  Valid  from  1  Aug  to 
30  Sep  2012. 


prescribed  observation  error,  the  /fit  value  should  be  at 
or  below  1. 

Figure  3  shows  the  /fit  value  for  the  temperature  (top 
panel),  salinity  (middle  panel),  and  velocity  (bottom 
panel)  observations  from  1  August  to  30  September, 
2012;  here  the  temperature  and  salinity  observations  are 
assimilated,  but  velocity  is  not.  Two  experiment  results 
are  shown  in  Fig.  3:  the  free-run  (FR)  NCOM  back¬ 
ground  (solid)  and  the  TS  analysis  (dash).  This  com¬ 
parison  is  made  to  determine  if  the  assimilation  of 
temperature  and  salinity  observations  improve  the 
analysis  compared  to  the  free-run  solution. 

Figure  3  shows  that  the  assimilation  of  temperature 
and  salinity  has  greatly  improved  the  observation  fit 
over  that  of  the  FR  background,  fitting  both  observation 
types  generally  within  the  prescribed  observation  error. 
The  fit  to  velocity,  here  an  independent  data  as  the  TS 
experiment  does  not  assimilate  velocity,  shows  no  im¬ 
provement  in  the  representation  of  the  surface  velocity 
over  that  of  the  FR  background.  As  one  would  expect, 
this  indicates  that  the  assimilation  of  temperature  and 
salinity  alone  is  not  enough  to  correct  the  ocean  velocity 
field.  One  can  examine  the  analysis  in  regards  to  with¬ 
held  observations  as  well.  Because  of  the  coarse  reso¬ 
lution  of  the  forecast  model  grid,  it  is  possible  to 
withhold  a  set  of  observations  for  each  data  type;  in  this 
case,  the  available  observations  are  thinned  by  removing 
those  observations  that  are  considered  redundant  (i.e., 


within  or  at  one  correlation  scale  distance  to  an  assimi¬ 
lated  observation).  A  subset  of  observations  can  be  col¬ 
lected  from  these  withheld  data  to  be  used  as  further 
validation  of  the  analysis.  Flere,  only  those  withheld 
observations  at  a  distance  of  at  least  30-40  km  from  an 
assimilated  observation  are  considered;  all  withheld 
observations  that  are  at  a  distance  less  than  this  threshold 
are  not  used  in  this  evaluation.  These  data  may  not 
necessarily  represent  purely  independent  observations; 
however,  evaluating  the  analysis  against  these  observa¬ 
tions  does  provide  insight  as  to  the  performance  of  the 
assimilation  far  from  the  location  of  included  observa¬ 
tions.  Figure  4  shows  the  /fit  value  of  FR  (solid)  back¬ 
ground  and  TS  (dash)  analysis  to  these  withheld 
observations.  The  fit  to  temperature  is  not  quite  as 
good  as  to  the  assimilated  observations,  however,  the 
TS  analysis  still  shows  some  improvement  over  the  FR 
background.  The  fit  to  salinity,  on  the  other  hand,  shows 
very  little  improvement  over  the  FR  background.  The 
error  in  salinity,  especially  during  September,  is  large 
with  /fit  values  exceeding  6.  A  further  check  of  these 
observations  shows  that  a  majority  of  them  are  located 
near  the  northern  Gulf  Coast,  in  the  Mississippi  River 
outflow  region;  these  observations  near  the  northern  Gulf 
Coast  are  more  numerous  in  the  month  of  September. 
These  specific  salinity  observations  exhibit  large  in¬ 
novations  when  compared  to  the  free-run  solution 
(some  greater  than  3.0  psu).  This  is  likely  due  to  poorly 
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Fig.  4.  As  in  Fig.  3,  but  for  nonassimilated  observations. 


specified  river  inflow  information  in  the  model  com¬ 
bined  with  incorrect  near-surface  ocean  currents,  as  seen 
in  the  velocity  /fit  values  in  Fig.  4.  It  is  important  to  note 
that  the  /fit  values  are  normalized  by  the  assumed  ob¬ 
servation  error,  which  in  this  case  is  0.1  psu  for  salinity. 
This  factor,  combined  with  the  large  innovations,  results 
in  a  large  /fit  value  for  salinity. 

Examining  the  forecasts  generated  from  the  4DVAR 
analysis  in  the  TS  experiment  can  be  done  by  computing 
the  skill  score  relative  to  the  FR  forecast.  The  skill  score 
is  a  measure  of  the  relative  root-mean-square  (RMS) 
error  in  one  forecast  solution  to  another.  Here,  the  RMS 
error  for  the  TS  and  FR  forecasts  are  computed  by  using 
all  of  the  available  temperature,  salinity,  and  velocity 
observations  valid  during  the  forecast  period: 


RMS  = 


M 

£  (y„ 


ra= 1 


Hjf 


(6) 


where  M  is  the  total  number  of  observations,  ym  is  the 
rath  observation,  and  is  the  model  forecast  solution.  It 
should  be  noted  that  the  observations  used  to  compute 
the  RMS  error  have  not  yet  been  assimilated  as  they  are 
compared  here  to  the  forecast,  not  the  analysis;  in  this 
regard,  the  data  used  in  this  comparison  are  independent 
observations.  The  skill  score  is  then  computed  by  com¬ 
paring  the  RMS  error  of  the  TS  and  FR  forecasts  as 


SS  =  1.0- 


RMSts 

RMSfr* 


(7) 


Equation  (7)  shows  that  if  the  RMS  error  of  the  TS 
forecast  is  lower  (higher)  than  the  FR  forecast,  the  skill 
score  metric  will  be  positive  (negative).  If  there  is  no 
change,  the  skill  score  value  should  be  nearly  zero. 
Figure  5  shows  the  computed  skill  score  from  the  TS 
experiment,  relative  to  the  FR  solution,  for  each  of  the 
96-h  forecasts  generated  from  the  4DVAR  analyses  ev¬ 
ery  4  days  for  temperature  (top  panel),  salinity  (middle 
panel),  and  velocity  (bottom  panel)  from  1  August  to 
30  September  2012  (solid  line  shows  the  skill  score  value, 
dashed  line  indicates  zero  skill  score  value).  Figure  5 
shows  that  the  TS  solution  produces  an  improved  fore¬ 
cast  in  temperature  and  salinity  during  most  of  the  ex¬ 
periment  time  frame  over  that  of  the  FR  solution.  For 
salinity,  the  TS  experiment  generally  outperforms  FR  in 
September,  with  a  few  exceptions.  There  is  almost  no 
improvement  whatsoever  in  the  velocity  forecast.  This  is 
not  surprising  given  the  analysis  results  shown  in  Figs.  3 
and  4.  Figure  5  suggests  that  the  assimilation  of  tem¬ 
perature  and  salinity  alone  does  not  constrain  the  model 
solution  enough  in  terms  of  the  surface  velocities. 

Adding  the  GLAD -derived  drifter  velocity  observa¬ 
tions  to  the  assimilation  is  done  in  the  ALL  experiment. 
One  concern  that  exists  when  combining  multiple  ob¬ 
servation  types  in  the  same  4DVAR  minimization  is  that 
one  or  more  data  types  will  be  fit  by  the  analysis  at  the 
expense  of  another  data  type.  There  are  many  reasons 
why  this  could  occur,  but  is  usually  due  to  the  relative 
errors  of  each  data  type,  the  size  of  the  innovations  and 
their  contribution  to  the  cost  function,  and  the  criteria 
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Fig.  5.  TS  forecast  skill  score  values,  measured  against  NCOM  free-run  solution  for  (top) 
temperature,  (middle)  salinity,  and  (bottom)  velocity.  Valid  from  1  Aug  to  30  Sep  2012.  Skill 
score  indicated  by  the  solid  line  and  zero  skill  score  value  indicated  by  the  dashed  line. 
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for  the  convergence  of  the  conjugate  gradient  algorithm. 
Temperature  observations,  in  particular,  can  be  numer¬ 
ous  (especially  at  the  surface)  in  comparison  to  other 
data  types  and  can  often  have  larger  innovations  than 
salinity  or  ocean  velocities.  Therefore,  if  the  convergence 


criterion  is  less  strict,  the  algorithm  may  converge  hav¬ 
ing  only  fitted  temperature  observations  properly.  For 
this  reason,  the  analysis  fit  to  temperature,  salinity,  and 
velocity  is  compared  to  the  TS  analysis  in  Fig.  6.  This 
figure  shows  the  /fit  value  to  determine  the  impact  of  the 
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Fig.  6.  The  /fit  metric  values  for  the  TS  analysis  solution  (solid  line)  and  ALL  analysis  so¬ 
lution  (dashed  line)  measured  against  (top)  assimilated  temperature  observations,  (middle) 
salinity  observations,  and  (bottom)  GLAD  velocity  observations.  Valid  from  1  Aug  to  30  Sep 
2012. 
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observations  on  this  analysis  and  compares  the  analysis 
error  in  temperature  (top  panel),  salinity  (middle  panel), 
and  velocity  (bottom  panel)  from  the  TS  experiment 
(solid)  and  the  ALL  experiment  (dashed).  Here,  it  is  seen 
that  the  ALL  experiment  continues  to  fit  the  temperature 
and  salinity  observations  well,  nearly  matching  the  results 
from  the  TS  experiment  (doing  marginally  better,  in  fact, 
than  the  TS  analysis  in  terms  of  salinity).  This  indicates 
that  the  inclusion  of  the  surface  velocity  data  type  in  the 
analysis  has  not  degraded  the  analysis  fit  to  the  other 
assimilated  data  types. 

Figure  6  also  shows  that  the  assimilation  of  surface 
velocity  from  the  GLAD  drifter  data  has  greatly  im¬ 
proved  the  analysis  fit  to  these  observations  over  the  TS 
analysis.  It  can  be  seen  that  the  ALL  experiment  fits  the 
velocity  observations  to  within  the  prescribed  error 
through  the  last  month  of  the  experiment,  although  it 
does  have  some  trouble  with  the  fit  earlier  through  the 
month  of  August.  This  is  likely  due  to  the  model  analysis 
and  forecast  adjusting  to  the  surface  currents  during  the 
earlier  portion  of  the  experiment. 

As  with  the  comparison  of  the  TS  analysis  and  FR 
background,  the  /fit  values  for  the  TS  and  ALL  analyses 
using  the  withheld  observations  are  computed;  Fig.  7 
shows  this  comparison.  The  /fit  value  for  temperature 
indicates  that  the  ALL  analysis  (dashed)  fits  the  with¬ 
held  temperature  observations  just  as  well  as  the  TS 
analysis  (solid).  In  addition  to  this,  the  ALL  analysis 
shows  an  improved  fit  to  the  withheld  salinity  observa¬ 
tions  when  compared  to  the  TS  analysis.  This  suggests 


that  the  assimilation  of  GLAD  velocity  observations 
near  the  northern  Gulf  Coast  is  helping  to  constrain  the 
salinity  gradients  near  the  Mississippi  River  outflow 
region,  which  results  in  a  better  fit  to  these  independent 
observations.  The  ALL  analysis  fit  to  withheld  GLAD 
velocity  observations  also  shows  improvement  over  the 
TS  analysis,  even  fitting  the  withheld  observations  to 
within  the  prescribed  data  error  during  the  latter  half  of 
the  experiment  run. 

Figure  8  examines  the  observation  minus  analysis 
misfit,  spatially,  and  by  velocity  component  from  1  Au¬ 
gust  to  30  September  for  the  u  component  (left  panels) 
and  v  component  (right  panels)  from  the  TS  analysis 
(top  panels)  and  ALL  analysis  (bottom  panels);  obser¬ 
vation  minus  analysis  misfits  are  plotted  daily  for  ease  of 
viewing.  Clearly,  the  ALL  analysis  fits  the  GLAD  ve¬ 
locity  data  better  than  the  TS  experiment,  with  most  of 
the  misfit  values  near  zero.  The  TS  analysis  appears  to 
have  higher  error  in  the  central  Gulf  of  Mexico  in  the 
vicinity  of  the  loop  current  eddy.  The  ALL  analysis 
appears  to  have  reduced  this  error  significantly  and  has 
a  better  representation  of  the  velocity  field  in  the  central 
Gulf  of  Mexico  than  the  TS  analysis. 

It  has  now  been  shown  that  the  inclusion  of  velocity 
observations  leads  to  an  improved  surface  velocity  rep¬ 
resentation  (according  to  available  observations)  than 
without  velocity  assimilation.  To  see  the  impact  of  these 
observations  on  the  subsequent  forecast,  Fig.  9  displays 
the  forecast  skill  score  of  the  ALL  96-h  forecasts,  relative 
to  the  TS  forecasts  using  (7)  for  temperature  (top  panel), 
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Fig.  8.  GLAD  velocity  observation  (top)  minus  TS  analysis  solution  and  (bottom)  minus  ALL  analysis  solution  for 
the  (left)  u  component  and  (right)  v  component  of  velocity.  Difference  values  plotted  daily  from  1  Aug  to  30  Sep 
2012.  Difference  values  in  ms-1  (magnitude  indicated  by  the  color  bar). 


salinity  (middle  panel),  and  velocity  (bottom  panel). 
The  ALL  forecast  of  temperature  is  just  as  good  as  the 
TS  solution,  indicating  that  the  ALL  forecast  is  also 
producing  an  improved  temperature  field  over  the  FR 


solution  (as  in  the  TS  experiment).  Interestingly,  the 
ALL  forecast  of  salinity  shows  a  significant  improve¬ 
ment  over  the  TS  forecast,  especially  through  the  month 
of  September.  This  is  consistent  with  the  results  in  Fig.  7 


Temperature  Skill  Score 


0.5 


Salinity  Skill  Score 


10  20  30  40  50  60 


Days  Past  August  L  2012 


Fig.  9.  ALL  forecast  skill  score  values,  measured  against  TS  forecast  solution  for  (top) 
temperature,  (middle)  salinity,  and  (bottom)  velocity.  Valid  from  1  Aug  to  30  Sep  2012.  Skill 
score  indicated  by  the  solid  line  and  zero  skill  score  value  indicated  by  the  dashed  line. 
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Fig.  10.  Observed  salinity  values  minus  (right)  TS  forecast  solution  and  (left)  ALL  forecast  solution,  valid  between 
10  and  30  Sep  2012  (values  in  psu,  indicated  by  color  bar).  NCOM  model  forecast  fields  are  mapped  to  the  observed 
location  at  the  appropriate  time. 


that  showed  the  ALL  analysis  fit  the  withheld  salinity 
observations  better  than  the  TS  analysis.  This  result 
suggests  that  the  better  fit  to  withheld  salinity  observa¬ 
tions,  along  with  the  improved  surface  current  represen¬ 
tation,  in  the  ALL  analysis  has  helped  to  improve  the 
prediction  of  salinity  in  the  ALL  forecast  relative  to 
the  TS  experiment.  The  bottom  panel  in  Fig.  9  displays 
the  forecast  skill  score  for  the  velocity  field  from  the  ALL 
experiment  relative  to  the  TS  experiment.  The  skill  score 
is  well  above  zero  for  the  majority  of  the  experiment  time 
frame,  indicating  that  the  improved  velocity  analysis  in 
the  ALL  experiment  does  indeed  translate  to  an  im¬ 
proved  forecast  when  compared  to  the  TS  experiment. 
The  improvement  gained  by  the  assimilation  of  the 
GLAD  velocity  observations  generally  lasts  the  entirety 
of  each  96-h  forecast,  indicating  that  the  improvement  is 
not  short  lived  and  the  memory  of  the  information  gained 
from  the  assimilation  in  the  forecast  is  significant. 

Examining  the  salinity  forecast  in  TS  and  ALL  more 
closely,  the  majority  of  the  improvement  in  the  skill  score 
in  Fig.  9  is  from  near-surface  salinity  observations  near 
the  northern  Gulf  Coast.  Figure  10  shows  the  absolute 


difference  between  a  set  of  near-surface  salinity  obser¬ 
vations  and  the  TS  salinity  forecast  (left  panel)  and  ALL 
salinity  forecast  (right  panel)  for  the  period  of  10-30 
September  2012  in  a  subset  region  from  the  NCOM  do¬ 
main  near  the  northern  Gulf  Coast  (absolute  difference 
value  indicated  by  the  color  bar).  It  is  clear  from  Fig.  10 
that  the  ALL  salinity  forecast  is  matching  the  observa¬ 
tions  much  more  closely  than  the  TS  forecast.  Further 
examining  the  forecast  salinity  within  the  same  region  as 
Fig.  10  (with  surface  velocities  overlaid)  from  25  Sep¬ 
tember  2012  from  TS  (Fig.  11,  left  panel)  and  ALL  (Fig. 
11,  right  panel)  shows  significant  differences  in  the  two 
forecast  salinity  fields.  The  salinity  forecast  from  TS 
shows  higher  salinity  values  south  of  28°N  with  the 
freshest  water  associated  with  the  Mississippi  River  out¬ 
flow  confined  nearest  the  coast.  This  differs  from  the 
ALL  salinity  forecast,  which  shows  generally  lower  sa¬ 
linity  values  south  of  28°N,  with  an  offshore  eddy  (cen¬ 
tered  at  26.5°N,  91°W)  driving  a  current  that  is  pulling 
freshwater  from  the  Mississippi  River  outflow  farther 
away  from  the  coast  into  the  central  Gulf.  It  was  shown  in 
Fig.  10  that  the  ALL  surface  currents  are  more  accurate 
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Fig.  11.  (left)  TS  forecast  salinity  field  and  (right)  ALL  forecast  salinity  field,  with  surface  velocity  field  overlaid 

(vector  plots),  valid  25  Sep  2012. 
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Fig.  12.  Absolute  dynamic  height  (ADH)  from  the  AVISO  product  (this  altimeter  product 
was  produced  by  Ssalto /Duacs  and  distributed  by  AVISO,  with  support  from  CNES  at  http:// 
www.aviso.oceanobs.com/duacs/),  valid  22  Aug  2012. 


than  TS  (when  compared  to  the  GLAD  velocity  obser¬ 
vations),  and  this  region  is  covered  quite  well  by  the 
GLAD  drifters  during  this  time  frame,  suggesting  that  the 
surface  currents  seen  in  the  ALL  forecast  in  Fig.  11  are 
likely  more  accurate  than  those  of  the  TS  forecast. 
Therefore,  it  is  likely  that  the  improved  near-surface 
velocity  forecast  in  ALL  is  driving  the  improved  salinity 
forecast  when  compared  to  the  TS  solution. 

The  assimilation  of  velocity  observations  can  also 
help  improve  the  mesoscale  eddy  representation  in  the 
model.  This  is  due  to  the  dynamical  balance  relationship 
provided  by  the  tangent  linear  and  adjoint  of  the  ocean 


model;  corrections  to  the  surface  velocity  field  can  lead  to 
an  improvement  in  the  model  surface  elevation  field. 
Figures  12,  13,  and  14  show  the  comparison  of  absolute 
dynamic  height  (ADH;  in  meters)  from  the  Archiving, 
Validation,  and  Interpretation  of  Satellite  Oceano¬ 
graphic  data  (AVISO)  product  [this  altimeter  product 
was  produced  by  the  Segment  Sol  multimissions 
d’Altimetrie,  d’ Orbit ographie  et  de  localisation  precise 
(Ssalto )/Data  Unification  and  Altimeter  Combination 
System  {Duacs)  and  distributed  by  AVISO,  with  support 
from  the  Centre  National  d’Etudes  Spatiales  (CNES)  at 
http://www.aviso.oceanobs.com/duacs/;  Fig.  12],  the  TS 
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Fig.  13.  As  in  Fig.  12,  but  for  absolute  dynamic  height  from  the  TS  forecast  solution. 
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Fig.  14.  As  in  Fig.  12,  but  for  absolute  dynamic  height  from  the  ALL  forecast  solution. 


experiment  (Fig.  13),  and  the  ALL  experiment  (Fig.  14), 
valid  on  22  August  2012.  When  comparing  TS  and  ALL 
to  the  AVISO  product  in  the  region  bounded  by  22°- 
30°N,  84°-92°W  (i.e.,  the  region  covered  by  dense  GLAD 
velocity  observations),  there  is  a  stark  contrast  between 
the  two  model  solutions.  The  TS  ADH  field  exhibits  one 
large,  elliptic  eddy  feature  south  of  the  Mississippi  River 
delta.  This  differs  from  the  AVISO  product  that  indicates 
the  eddy  structure  is  smaller  with  an  elongated  “tongue” 
of  high  ADH  extending  to  the  south  and  east  of  the  pri¬ 
mary  eddy  (with  regions  of  lower  ADH  to  the  north  and 
south  of  this  tongue).  The  ADH  field  from  the  ALL 
forecast  exhibits  an  eddy  structure  much  closer  to  the 
AVISO  product  than  does  the  TS  forecast.  The  ALL 
ADH  field  shows  the  same  tongue  structure  with  the 
lower  ADH  regions  north  and  south  as  in  the  AVISO 
product,  albeit  with  higher  overall  height  values.  This 
indicates  that  the  assimilation  of  surface  velocity  mea¬ 
surements  from  GLAD  has  helped  to  constrain  the 
surface  eddy  field  in  the  vicinity  of  the  drifter  observa¬ 
tions.  An  examination  of  other  time  periods  throughout 
the  ALL  experiment  show  similar  agreement  with  the 
AVISO  product  (figures  not  shown). 

6.  Summary 

The  assimilation  of  surface  velocity  observations, 
derived  from  GLAD  drifter  positions,  within  the  weak- 
constraint  NCOM-4DVAR  has  been  shown.  By  com¬ 
paring  the  assimilation  results  in  terms  of  the  analysis 
fit  to  the  observations,  as  well  as  the  subsequent  forecast 
fit  to  future  observations  between  an  experiment  with 
no  velocity  assimilation  (TS)  and  an  experiment  with 


velocity  assimilation  (ALL),  we  have  been  able  to  con¬ 
clude  that  assimilating  the  surface  velocity  observations 
leads  to  a  substantial  improvement  in  not  only  the 
analysis  fit,  but  also  in  the  forecast  as  well.  The  assimi¬ 
lation  of  velocity  observations  also  led  to  an  improved 
salinity  forecast,  as  shown  in  the  comparison  of  the  re¬ 
sults  from  the  ALL  and  TS  experiments.  This  is  likely 
due  to  the  fact  that  the  assimilation  of  GLAD  velocity 
observations  near  the  northern  Gulf  Coast  is  helping  to 
constrain  the  salinity  gradients  near  the  Mississippi 
River  outflow  region,  resulting  in  a  better  fit  to  these 
salinity  observations  in  that  region.  The  forecast  skill 
score  of  the  ALL  forecast  relative  to  the  TS  forecast 
shows  that  the  information  gained  in  the  assimilation  of 
velocity  observations  is  not  short  lived,  with  improve¬ 
ment  in  the  velocity  forecast  in  ALL  over  the  TS  ex¬ 
periment  shown  to  remain  out  to  96  h.  It  was  also  shown 
that  the  assimilation  of  surface  velocity  measurements 
can  improve  the  model  representation  of  the  surface 
elevation  field.  A  comparison  to  the  AVISO  absolute 
dynamic  height  product  shows  that  the  ALL  experiment 
captures  the  mesoscale  eddy  structure  more  accurately 
than  the  TS  experiment  in  the  vicinity  of  the  velocity 
measurements.  In  short,  it  was  shown  that  the  assimila¬ 
tion  of  temperature  and  salinity  observations  alone  in 
the  4DVAR  does  not  properly  constrain  the  model 
representation  of  the  surface  flow  field  in  the  Gulf  of 
Mexico.  The  results  of  this  study  suggest  that  surface 
velocity  measurements  can  be  used  to  help  constrain  the 
model  flow  field  and,  therefore,  the  mesoscale  eddy 
structures.  Assimilating  surface  velocity  information 
can  help  to  correct  the  model  representation  of  the 
surface  elevation  field  as  well  as  help  to  constrain  the 
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tracer  fields,  specifically  the  surface  salinity.  The  results 
shown  in  this  study  suggest  that  surface  velocity  obser¬ 
vations  can  have  a  profound  impact  on  the  model  solu¬ 
tion,  not  only  on  the  representation  of  the  flow  field,  but 
across  other  ocean  variables  as  well. 

Now  that  the  assimilation  of  surface  velocity  obser¬ 
vations  has  been  shown  to  work  well  when  the  velocity 
observations  take  the  Eulerian  form,  attention  can  be 
given  to  the  Lagrangian  assimilation  of  the  GLAD  drifter 
observations  in  their  native  Lagrangian  form.  To  do  this, 
the  capability  to  evolve  model  drifters  within  the  NCOM 
solution  will  be  added  in  the  near  future.  Once  this  is 
done,  the  NCOM-4DVAR  can  be  used  to  minimize  the 
distance  between  the  model  and  observed  drifter  posi¬ 
tions  (using  a  suitable  measurement  functional),  thereby 
propagating  this  information  to  the  model  Eulerian  ve¬ 
locity  fields  in  the  analysis  step.  The  results  of  a  future 
Lagrangian  data  assimilation  effort  with  GLAD  data  will 
be  compared  to  this  work  to  determine  if  there  is  any 
improvement  in  the  accuracy  of  the  velocity  analysis  and/or 
forecast  over  the  pseudo-Lagrangian  DA  method,  and 
if  so,  it  will  determine  if  the  gain  is  substantial.  This 
question  will  be  investigated  in  a  future  experiment  and 
subsequent  paper. 
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